MEDB 5501, Module03

---
title: "Normal probabilities and quantiles"
format: 
  revealjs:
    slide-number: true
    embed-resources: true
editor: source
---

This program displays the standard normal curve along with various probabilities and quantiles. It was written by Steve Simon on 2024-09-01 and is placed in the public domain.

## Load the tidyverse library

```{r setup}
#| message: false
#| warning: false
library(tidyverse)
```

## Using R to draw the standard normal curve

use seq to calculate 100 evenly spaced values between -4 and +4 and dnorm to compute the bell curve at each point. Use geom_polygon to paint the area surrounded by the bell curve.

```{r standard-normal}
x <- seq(-4, 4, length=100)
y <- dnorm(x)
data.frame(x, y) |>
  ggplot(aes(x, y)) + 
    geom_polygon(fill="white", color="black") -> normal_curve
normal_curve
```

## P[Z < 1.5]

Use geom_vline and geom_label to draw a vertical reference line and add text to the normal curve. The pnorm function computes the standard normal probability.

```{r prob-1}
a <- 1.5
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label=a) +
  geom_label(x=a-0.5, y=0, label="Area = ?")
pnorm(1.5)
```

## P[Z < -0.5]

```{r prob-2}
a <- -0.5
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label=a) +
  geom_label(x=a-0.5, y=0, label="Area = ?")
pnorm(-0.5)
```

## P[Z > 1]

When you are calculating the probability on the right (probability greater than some number), use 1-pnorm.

```{r prob-3}
a <- 1
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label=a) +
  geom_label(x=a+0.5, y=0, label="Area = ?")
1- pnorm(1)
```

## P[Z > -2]

```{r prob-4}
a <- -2
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label=a) +
  geom_label(x=a+0.5, y=0, label="Area = ?")
1- pnorm(-2)
```

## P[-2.5 < Z < 2.5]

When you are calculating the probability between two values, compute pnorm of the larger value minus pnorm of the smaller value.

```{r prob-5a}
a <- 2.5
normal_curve +
  geom_vline(xintercept=-a) +
  geom_vline(xintercept= a) +
  geom_label(x=-a, y=0.4, label=-a) +
  geom_label(x= a, y=0.4, label= a) +
  geom_label(x=0, y=0, label="Area = ?")
pnorm(2.5)- pnorm(-2.5)
```

## P[-0.5 < Z < 0.5]

```{r prob-6}
a <- 0.5
normal_curve +
  geom_vline(xintercept=-a) +
  geom_vline(xintercept= a) +
  geom_label(x=-a, y=0.4, label=-a) +
  geom_label(x= a, y=0.4, label= a) +
  geom_label(x=0, y=0, label="Area = ?")
pnorm(0.5) - pnorm(-0.5)
```

## 25th percentile of a standard normal

Use qnorm to calculate quantiles of the standard normal distribution.

```{r quantile-1}
p <- 0.25
a <- qnorm(p)
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label="Quantile = ?") +
  geom_label(x=a-0.5, y=0, label=p)
qnorm(0.25)
```

## 90th percentile of a standard normal

```{r quantile-2}
p <- 0.9
a <- qnorm(p)
normal_curve +
  geom_vline(xintercept=a) +
  geom_label(x=a, y=0.4, label="Quantile = ?") +
  geom_label(x=a-0.5, y=0, label=p)
qnorm(0.9)
```

   x r
1  7 4
2  3 2
3 23 9
4  2 1
5  5 3
6 13 6
7 11 5
8 17 7
9 19 8

   x r pctile
1  7 4    0.4
2  3 2    0.2
3 23 9    0.9
4  2 1    0.1
5  5 3    0.3
6 13 6    0.6
7 11 5    0.5
8 17 7    0.7
9 19 8    0.8

   x r pctile          z
1  7 4    0.4 -0.2533471
2  3 2    0.2 -0.8416212
3 23 9    0.9  1.2815516
4  2 1    0.1 -1.2815516
5  5 3    0.3 -0.5244005
6 13 6    0.6  0.2533471
7 11 5    0.5  0.0000000
8 17 7    0.7  0.5244005
9 19 8    0.8  0.8416212

---
data_dictionary: fev (.csv, sas7bdat, .sav, .txt)

copyright: >
  The author of the jse article holds the copyright, but does not list
  conditions under which it can be used. Individual use for educational
  purposes is probably permitted under the Fair Use provisions of
  U.S. Copyright laws.

description: >
  Forced Expiratory Volume (FEV) in children. The data was collected 
  in Boston in the 1970s.

additional_description:
  https://jse.amstat.org/v13n2/datasets.kahn.html

download_url:
  https://www.amstat.org/publications/jse/datasets/fev.dat.txt

format:
  csv: comma delimited
  sas7bdat: proprietary (SAS)
  sav: proprietary (SPSS)
  txt: fixed width
  

varnames: 
  not included
  

missing_value_code:
  not needed
  

size:   
  rows: 654
  columns: 5

vars:
  age:
    scale: ratio
    range: positive integer
    unit: years
    

  fev:
    label: Forced Expiratory Volume
    scale: ratio
    range: positive real
    unit: liters
    

  ht:
    label: Height
    scale: positive real
    unit: inches

  sex:
    value:
      0: Female
      1: Male

  smoke:
    value:
      0: Nonsmoker
      1: Smoker
---

---
title: "Analysis of fev data"
format: 
  html:
    embed-resources: true
editor: source
---

This program assesses the normality of variables in a study of pulmonary function in children. There is a [data dictionary][dd] that provides more details about the data. The program was written by Steve Simon on 2024-09-02 and is placed in the public domain.

[dd]: https://github.com/pmean/datasets/blob/master/fev.yaml

## Libraries

The tidyverse library is the only one you need for  this program.

```{r setup}
#| message: false
#| warning: false
library(tidyverse)
```

## List variable names

Since the variable names are not listed in the data file itself, you need to list them here.

```{r names}
fev_names <- c(
    "age",
    "fev",
    "ht",
    "sex",
    "smoke")
```

## Reading the data

Here is the code to read the data and show a glimpse. 

```{r read}
fev <- read_csv(
  file="../data/fev.csv",
  col_names=fev_names,
  col_types="nnncc")
glimpse(fev)
```

## Calculate mean and standard deviation for fev

To orient yourself to the data, calculate a few descriptive statistics.

```{r descriptive-fev}
fev |>
  summarize(
    fev_mean=mean(fev),
    fev_stdv=sd(fev))
```

## Histogram for fev, wide bars

```{r histogram-fev-wide}
ggplot(data=fev, aes(x=fev)) +
  geom_histogram(
    binwidth=0.5,
    color="black",
    fill="white")
```

## Histogram for fev, narrow bars

```{r histogram-fev-narrow}
ggplot(data=fev, aes(x=fev)) +
  geom_histogram(
    binwidth=0.1,
    color="black",
    fill="white")
```

Although some may interpret these histograms as showing a slight skewness, I would interpret them as being approximately normal.

## Normal probability plot for fev

The qqnorm function produces a normal probability plot. The default option for most plots is landscape orientation (the width is larger than the height). The q-q plot, however, looks best if figure width and height are equal.

```{r qqplot-fev}
#| fig-width: 5
#| fig.height: 5
qqnorm(fev$fev)
```

The normal probability plot is reasonably close to a straight line, indicating that the data comes reasonably close to following a normal distribution.